Spatial modeling reveals nuclear phosphorylation and subcellular shuttling of YAP upon drug-induced liver injury

The Hippo signaling pathway controls cell proliferation and tissue regeneration via its transcriptional effectors yes-associated protein (YAP) and transcriptional coactivator with PDZ-binding motif (TAZ). The canonical pathway topology is characterized by sequential phosphorylation of kinases in the cytoplasm that defines the subcellular localization of YAP and TAZ. However, the molecular mechanisms controlling the nuclear/cytoplasmic shuttling dynamics of both factors under physiological and tissue-damaging conditions are poorly understood. By implementing experimental in vitro data, partial differential equation modeling, as well as automated image analysis, we demonstrate that nuclear phosphorylation contributes to differences between YAP and TAZ localization in the nucleus and cytoplasm. Treatment of hepatocyte-derived cells with hepatotoxic acetaminophen (APAP) induces a biphasic protein phosphorylation eventually leading to nuclear protein enrichment of YAP but not TAZ. APAP-dependent regulation of nuclear/cytoplasmic YAP shuttling is not an unspecific cellular response but relies on the sequential induction of reactive oxygen species (ROS), RAC-alpha serine/threonine-protein kinase (AKT, synonym: protein kinase B), as well as elevated nuclear interaction between YAP and AKT. Mouse experiments confirm this sequence of events illustrated by the expression of ROS-, AKT-, and YAP-specific gene signatures upon APAP administration. In summary, our data illustrate the importance of nuclear processes in the regulation of Hippo pathway activity. YAP and TAZ exhibit different shuttling dynamics, which explains distinct cellular responses of both factors under physiological and tissue-damaging conditions.


Introduction
The evolutionary conserved Hippo signaling pathway controls tissue homeostasis through the regula-

Results
Establishment of an in vitro model for measuring time-resolved spatial localization of YAP and TAZ in hepatocellular cells Due to technical limitations, primary human hepatocytes are not suitable for the analysis of dynamic YAP/TAZ shutting in vitro: they rapidly undergo trans-differentiation in culture and different human donors exhibit genomic/epigenetic variability, which diminish reproducibility. We therefore selected the hepatocyte-derived liver cancer cell line Hep3B that was characterized by a prominent YAP/TAZ and cytochrome P450 2E1 (CYP2E1) expression (Figure 1-figure supplement 1A; Weiler et al., 2020). Detectable CYP2E1 levels are required for the enzymatic transformation of APAP to toxic N-acetyl-p-benzoquinone imine (NAPQI) (Lee et al., 1996).
Since we were interested in cell-to-cell variability and spatial localization information of YAP and TAZ, we established a cell line that allowed the quantitative and spatiotemporal investigation of subcellular YAP/TAZ distribution near the single-cell resolution. For this, we generated Hep3B cells that stably express mVenus-tagged YAP (mVenus-YAP) and mCherry-tagged TAZ (mCherry-TAZ) ( Figure 1A). In addition, mCerulean-tagged histone H2B (mCerulean-H2B) allowed the spatial allocation of YAP/TAZ. These chosen fluorophores facilitated optimal discrimination between emission spectra by live cell microscopy ( Figure 1-figure supplement 1B). Subcellular localization changes of tagged YAP and TAZ under variable cell density conditions illustrated functionality of both proteins regarding nuclearcytoplasmic shuttling ( Figure 1B). However, we observed a high degree of intracellular variability.
For analysis of variable high-throughput data, we developed an algorithm for the semiquantitative measurement of tagged YAP and TAZ in nuclei and cytoplasm. In brief, spatial image data were classified based on the expression of mVenus-YAP and mCherry-TAZ using the Weka segmentation algorithm ( Figure 1C, Figure 1-figure supplement 1C; Arganda-Carreras et al., 2017). The fluorophore intensity measurements were used to calculate the nuclear/cytoplasmic ratio (NCR), which characterized the relative nuclear enrichment of YAP or TAZ. Therefore, this approach allowed us to simultaneously define slightest subcellular changes of YAP and TAZ in small cell populations (field of view, FOV) as well as individual cells under live cell conditions. Subsequently, the fluorophore-expressing cells were grown under variable cell density conditions and the YAP/TAZ localization was measured by confocal microscopy. The results confirmed that tagged YAP and TAZ dynamically shuttled between cytoplasm and the nucleus depending on abundance of cell-cell contacts (cell density) ( Figure 1D). However, TAZ showed a considerably less pronounced dynamic behavior compared to YAP ( Figure 1D). YAP was predominantly detectable in cell nuclei under low cell density conditions (NCR > 1), while increasing cell density led to cytoplasmic enrichment of both factors (NCR < 1). Using this algorithm, cell density-dependent protein shuttling was also demonstrated for other liver cancer cells (Tóth et al., 2021).
Together, the confocal imaging and image processing pipeline efficiently detects subtle differences in the dynamic subcellular changes of YAP and TAZ.

Mathematical modeling predicts nuclear phosphorylation of YAP and TAZ
To investigate why Hippo pathway effectors YAP and TAZ differently respond to a low cell density and which Hippo pathway topology can explain the observed localization differences, we mathematically modeled the Hippo pathway using a PDE modeling framework (for detailed description of the PDE modeling approach refer to Materials and methods and Appendix 1).
First, we investigated a canonical Hippo pathway model, where phosphorylation takes place in the cytoplasm and where unphosphorylated YAP/TAZ shuttle to the nucleus ( Figure 2A). However, this mathematical model was not able to reproduce the localization patterns which were observed in the experimental data ( Figure 2B). For example, the model cannot sufficiently explain the distribution of YAP and TAZ in cell nuclei or at the nuclear membrane. In detail, the cell area around the nuclear envelope on the cytoplasmic side strongly underrepresented the concentration of YAP/TAZ, as indicated by the blue pixels in the residual image (experimental data minus model simulation). Moreover, the simulated YAP/TAZ distribution pattern of the canonical model showed a decreased protein concentration in the center of the nucleus, with increasing gradient toward the nuclear envelope. This indicates that the experimental observation is dominated by a protein that is disseminated from the The hepatocyte-derived cell line Hep3B was transduced using three lentiviral vectors coding for mCerulean-tagged H2B (pLentiPGK-mCerulean-H2B), mVenus-tagged YAP (pRRLN-EF1α-mVenus-YAP), and mCherry-tagged TAZ (pRRLN-EF1α-d2mCherry-TAZ). Combined treatment with the antibiotics hygromycin, geneticin, and blasticidin selected for triple-positive cells. Cells were confocally imaged in three channels: 445 nm (CFP, mCerulean), 514 nm (YFP, mVenus), and 561 nm (RFP, mCherry). Scale bar: 250 µm. (B) Exemplary pictures illustrating the subcellular localization of H2B, YAP, and TAZ proteins under low (left) and high (right) cell density conditions. Scale bar: 50 µm. (C) Automatic image analysis workflow depicts analysis of nuclear (left) and cytoplasmic (right) fluorescence intensity in confocal images of living cells. The acquired images were prescreened for imaging artifacts (e.g., outof-focus images were excluded) and subjected to the image analysis pipeline in Fiji. Object classification (Weka algorithm), thresholding, and object counting were performed. Object masks were overlaid with original data and the nuclear/cytoplasmic ratio (NCR) was calculated by dividing average fluorescent signal intensity of YAP/TAZ in the nucleus with the fluorescence intensity in the cytoplasm. Scale bar: 50 µm. (D) Quantification of YAP and TAZ NCR under increasing cell density conditions. One of four representative experiment is depicted. Left: the NCR of nucleus and undergoes diffusion and degradation. The canonical Hippo pathway cannot explain this effect, illustrated by residuals ( Figure 2B). This let us conclude that alternative model topologies must be considered to explain the experimental data. We therefore generated an alternative PDE model, which precisely described YAP/TAZ localization as observed in the image data. In this model, a phosphorylation and dephosphorylation reaction of YAP and TAZ in the nucleus was included ( Figure 2C,D, Figure 2-figure supplement  1A, B). Moreover, the alternative model describes that unphosphorylated YAP/TAZ are transported to the nucleus and pYAP/pTAZ are rapidly excluded from the nucleus, which agrees with the general consensus that phosphorylated YAP and TAZ predominantly localize in the cytoplasm (Figure 2figure supplement 1C). Using well-fitting parametrization of the model, the simulated concentration pattern was dominated by phosphorylated YAP/TAZ disseminating from the nucleus, which now matched the observed experimental data well ( Figure 2D).
Since YAP and TAZ are considered biochemically similar molecules, the question arises how the pronounced observed differences in the spatial distribution patterns can be explained. Using the alternative Hippo pathway model, we demonstrated that the different phenotypes could be reproduced with one single PDE model where only the parameters of the YAP/TAZ phosphorylation and dephosphorylation reactions in the nucleus differed, while all other parameters (including nuclear import and export reactions) were kept the same ( Figure 2E, Figure 2-figure supplement 1D).
In detail, low nuclear phosphorylation/dephosphorylation ratio (0.17) induced nuclear localization of YAP, while high ratio (0.75) induced nuclear protein exclusion. Upon increased phosphorylation activity (>0.2), the nuclear fraction of TAZ decreased and the summarized residual term for TAZ was reduced. According to residual terms, best accordance between the model and the experimental data was achieved at low phosphorylation/dephosphorylation ratio for YAP and at higher phosphorylation/ dephosphorylation ratio for TAZ ( Figure 2E, right panel). These findings suggest that differences in the nuclear phosphorylation rates of YAP and TAZ can explain the distinct shuttling patterns of both factors (as observed in Figure 1B). In addition, TAZ may require stronger phosphorylation than YAP to induce its efficient shuttling.
The results of the alternative mathematical model required that major components for YAP/TAZ phosphorylation process were present in the cell nucleus. Indeed, subcellular fractionation experiments illustrated that the Hippo pathway kinases LATS1/2 (total and phosphorylated) were predominantly localized in nuclear protein fraction ( Figure 2F). Moreover, we and others reported that nuclear LATS1/2 can control YAP phosphorylation and localization (data not shown) (Ege et al., 2018;Li et al., 2014). The relevance of this nuclear interaction is substantiated by a proximity ligation assay (PLA), which illustrated that pLATS1/2 physically interacted with YAP not only in the cytoplasm but also in the nucleus ( Figure 2G). To quantitatively compare the number of YAP/pLATS interactions in cell nuclei and cytoplasm, we utilized a computational algorithm counting signals in nuclear and cytoplasmic area of the cells (Figure 2-figure supplement 2A, B). The signal quantification confirmed that YAP and pLATS interactions were detectable in both compartments with significantly higher interaction count in cell nuclei ( Figure 2H).
In summary, our findings suggest that nuclear YAP/TAZ phosphorylation contributes to their inactivation and nuclear exclusion. Moreover, variable phosphorylation rates of YAP and TAZ explain the subcellular shuttling differences between both factors.

APAP regulates YAP protein localization and activity
To test if our findings are of relevance under conditions where YAP and TAZ are actively regulated, we decided to use a drug-induced liver injury (DILI) model. For this, we treated Hep3B cells expressing YAP (yellow, n = 310) and TAZ (red, n = 310) was plotted against cell count per visual field (0.4 mm 2 ). A high NCR value indicates predominant nuclear localization. Dashed line shows mean cell count. Right: violin plots summarize shift in NCR between low (below mean cell count) and high (above mean cell count) cell density for YAP and TAZ.
The online version of this article includes the following source data and figure supplement(s) for figure 1:    Image quantification revealed that APAP led to a gradual and time-dependent nuclear enrichment of YAP after 48 hr, while TAZ only weakly responded to APAP ( Figure 3A,B). This APAP-induced nuclear shuttling effect was clearly detectable under high cell density culture conditions, which is characterized by nuclear YAP/TAZ exclusion. Effects on YAP were less pronounced 24 hr after APAP administration ( Figure 3-figure supplement 1A ,B). No obvious response was detectable for YAP and TAZ at earlier time points post APAP treatment (data not shown). Severe effects of APAP on cell toxicity and apoptosis in the chosen experimental setup were excluded by measuring cell viability and PARP cleavage (data not shown).
The APAP-induced nuclear enrichment of YAP was confirmed with cell population-based Western immunoblotting followed by detection of total and phosphorylated YAP ( Figure 3C, for quantification of all western blots see Figure 3-figure supplement 2). Notably, at early time points YAP hyperphosphorylation indicated protein inactivation (up to 3 hr after APAP administration). However, between 24 and 48 hr after APAP treatment, a clear YAP dephosphorylation was detectable ( Figure 3C). Similar results were observed for another hepatocyte-derived cell line ( Figure 3-figure supplement 1C).
To test if APAP-dependent dephosphorylation/activation of YAP at later time points is leading to its transcriptional activation, we measured the relative expression of known YAP target genes cysteine-rich angiogenic inducer 61 (CYR61) and ankyrin repeat domain 1 (ANKRD1) by real-time PCR. Indeed, we observed induction of these genes already 24 hr after APAP treatment, which reflected the earliest activation of YAP ( Figure 3D). To confirm YAP dependence on the target gene induction and to exclude YAP-independent effects on gene expression, we performed additional rescue experiments. For this, the expression of YAP was silenced by RNA interference (RNAi) in APAP-treated cells followed by the measurement of YAP target genes. According to our hypothesis, YAP inhibition partly abolished the APAP-dependent induction of CYR61 and ANKRD1 ( Figure 3E).
Together, these data demonstrate that APAP predominantly acts on the Hippo pathway effector YAP in a bimodal manner. Upon APAP treatment an immediate phosphorylation/inactivation of YAP is followed by its late dephosphorylation/activation.

APAP controls YAP phosphorylation via ROS and AKT
Based on previous data, we hypothesized a mechanistic connection between APAP-induced ROS (Barbier-Torres et al., 2017;Shuhendler et al., 2014) and AKT-driven YAP phosphorylation (Basu et al., 2003;Romano et al., 2010). This mechanistic link was investigated at early time points after APAP treatment (up to 3 hr) to exclude unspecific effects caused by APAP at later time points (e.g., due to secondary and/or unspecific APAP effects). rapidly exported to the cytoplasm. (D) PDE simulation of the alternative model compared to experimental data. Residuals for the subcellular localization (e.g., nuclear distribution) sufficiently reflect results from confocal microscopy. Nuclear outline is indicated in black. (E) Simulated impact of the nuclear phosphorylation to dephosphorylation ratio on subcellular localization of YAP in the alternative Hippo pathway model. Left: model simulation of two phosphorylation/dephosphorylation rates (0.17 and 0.75). Right: summarized residuals with respect to experimental data of YAP and TAZ as a function of phosphorylation to dephosphorylation ratio in the nucleus. (F) Western immunoblot after nuclear and cytoplasmic protein fractionation under low cell density conditions. The central Hippo pathway kinases LATS1/2 were detectable in the nuclear fraction. As expected, high pYAP levels are detectable in the cytoplasm, illustrating that the protein is transported outside the nucleus upon phosphorylation. PARP and Tubulin serve as loading controls for nuclear and cytoplasmic fractions, respectively. Equal amounts of cytoplasmic and nuclear proteins were loaded (n = 4; one representative experiment shown). (G) Proximity ligation assay (PLA) for phosphorylated LATS1/2 (pLATS1/2) and YAP (top row). Red dots indicate physical interaction between pLATS1/2 and YAP. DAPI (4′,6-diamidino-2-phenylindole)-stained nuclei are indicated in cyan. Individual pLATS and YAP antibodies serve as assay controls. Bottom row: computational segmentation of nuclei (empty circles) and dots (black spots). Scale bar: 50 µm. (H) Quantification of the PLA assay. Each dot represents an individual image. Left: quantification of interactions between YAP and pLATS1/2 in the nucleus (n = 18) compared to the negative controls (pLATS1/2 and YAP antibodies alone, n = 16 and n = 17, respectively). Right: the number of PLA dots of pLATS1/2 and YAP interaction in cytoplasm and in nuclei (n = 18). Statistical test: two-tailed paired t-test (p value = 0.01) **p ≤ 0.01.
The online version of this article includes the following source data and figure supplement(s) for figure 2: Source data 1. Raw data for Western Blot shown in Figure 2F.   Western immunoblotting showed that both YAP and AKT were phosphorylated early after APAP incubation ( Figure 4A), which was in accordance with our previous findings ( Figure 3C). In parallel, an ROS detection assay revealed that APAP-induced prominent ROS activity comparable to the known ROS inducers hydrogen peroxide (H 2 O 2 ) and tert-butyl hydroperoxide (TBHP) in hepatocellular Figure 3. APAP regulates YAP phosphorylation, localization, and expression of its target genes. (A) Live cell confocal imaging of H2B-mCerulean, YAP-mVenus, and TAZ-mCherry in Hep3B cells. Upper row: control treatment (phosphate-buffered saline, PBS). Lower row: APAP treatment (10 mM) induces nuclear enrichment of YAP, but not TAZ protein after APAP treatment within 48 hr (n = 4; one representative experiment shown). Scale bar: 50 µm. (B) Live cell confocal microscopy image quantification. Nuclear/cytoplasmic ratio (NCR) of YAP and TAZ with APAP (10 mM for 48 hr; yellow and red circles, n = 180 per channel) and control (PBS; gray circles, n = 180 per channel) treatment under increasing cell density conditions (represented as mean cell count per visual field). Each dot represents a single visual field. Dashed line: mean cell density. (C) Western immunoblot of YAP and pYAP after APAP treatment (10 mM) in Hep3B cells (n = 6; one representative experiment shown). GAPDH served as a loading control. (D) Relative expression of the YAP target genes ANKRD1 and CYR61 24 hr after APAP (10 mM) treatment in HLF cells (n = 2; one out of two biological replicates shown). (E) Rescue experiment: relative expression of YAP and its target genes CYR61 and ANKRD1 with and without siRNA-mediated YAP inhibition (for 24 hr) with or without APAP (10 mM) treatment for 24 hr (n = 2; one out of two biological replicates shown).
The online version of this article includes the following source data and figure supplement(s) for figure 3: Source data 1. Source data for Western Blot shown in Figure 3C.    . These data indicated that APAP regulated ROS activity, which itself controlled AKT and YAP activity.
To confirm the rapid and AKT-dependent phosphorylation of YAP, liver cells were treated with hepatocyte growth factor (HGF), which activates/phosphorylates AKT kinases via binding the receptor c-MET within 10 min (Xiao et al., 2001). HGF administration led to a clear phosphorylation of AKT but not YAP, which might be caused by saturation effects (no further YAP phosphorylation possible under given culture conditions) ( Figure 4E). However, simultaneous inhibition of AKT kinase activity by the specific inhibitor AKTi (Adlung et al., 2017;Barnett et al., 2005), not only abolished AKT phosphorylation but also reduced YAP phosphorylation ( Figure 4E). As expected, upon AKTi administration concomitantly with APAP, the phosphorylation of both AKT and YAP decreased ( Figure 4F).
Since our alternative PDE model predicted that the nuclear YAP phosphorylation is crucial for protein shuttling and its activity ( Figure 2C-E), we hypothesized that AKT and YAP not only physically interact in the cytoplasm but also in cell nuclei. Indeed, PLA illustrated that AKT and YAP, as well as the respective phosphorylated isoforms, interacted in the nucleus ( Figure 4G,H). Although, interaction between AKT/YAP and pAKT/pYAP was also detectable in the cytoplasm, a prominent protein colocalization took place in the nuclear compartment ( Figure 4G). Importantly, image-based quantification revealed that the interaction between nuclear YAP and AKT increased upon 6hr treatment with APAP ( Figure 4I), indicating an APAP-dependent shift to nuclear phosphorylation and interaction, as predicted by the PDE model. Indeed, blocking AKT activity by AKTi reduced YAP phosphorylation, moderately elevated nuclear YAP positivity, and increased the physical interaction between AKT and YAP in cell nuclei ( Figure 4F and data not shown).
In summary, our data show that APAP affects YAP activity and shuttling behavior by enhancing cellular ROS followed by nuclear AKT/YAP binding and YAP phosphorylation.

Sequential activation of ROS, AKT, and Hippo/YAP in mouse livers after APAP intoxication
The results of our mathematical model ( Figure 2) and the in vitro experiments (Figure 3 and Figure 4) demonstrated a distinct sequence of molecular events that control YAP activity after APAP exposure: ROS induction is followed by AKT activation and YAP phosphorylation/inactivation (early events -up to 6 hr). This is followed by phase, where YAP is dephosphorylated and transcriptionally active (late events -24 to 48 hr).
To confirm this in vivo, mice were injected with a hepatotoxic dose of APAP (300 mg/kg) and liver tissues were collected up to 16 days (in total at nine time points). Subsequently, liver specimens of two biological replicates is shown). APAP as well as H 2 O 2 and TBHP induce ROS formation in living cells. Statistical test: one-way analysis of variance (ANOVA) with Geisser-Greenhouse correction (adjusted p value = 0.001) ***p-value ≤ 0.001. Whiskers depict min and max of the dataset. (C) Western immunoblot analysis of pYAP, YAP, pAKT, and AKT after H 2 O 2 (2 mM) and APAP (10 mM) treatment for 1 hr. APAP, as well as H 2 O 2 induce YAP and AKT phosphorylation compared to untreated cells (ctrl); n = 4. (D) Western immunoblot analysis of pYAP, YAP, pAKT, and AKT dynamics after TBHP (300 µM) and APAP (10 mM) treatment for 2 hr (Hep3B). APAP and TBHP induce AKT and YAP phosphorylation as compared to respective controls (ctrl). (E) Western immunoblot analysis of pAKT, AKT, pYAP, and YAP after hepatocyte growth factor (HGF) treatment (10 ng/µl) for 10 min (Hep3B). The cells were pretreated with an AKT inhibitor VIII (AKTi, 1, 5, and 10 µM) and starved for 3 hr prior HGF administration. Data illustrate that AKT inhibition prevents HGF-induced AKT and YAP phosphorylation. (F) Western immunoblot of pAKT, AKT, pYAP, and YAP after APAP (10 mM) administration for 1 hr (Hep3B). Cells were starved in Fetal calf serum (FCS)-free medium and pretreated with AKTi (5 µM) before APAP treatment for 3 hr. Results demonstrate that AKT inhibition reduces YAP phosphorylation early after APAP treatment. (G) Proximity ligation assay (PLA) of the protein combinations YAP/AKT and pYAP/ pAKT. Red dots indicate interactions between YAP/AKT and pYAP/pAKT, respectively. DAPI-stained nuclei are depicted in cyan. Single YAP, AKT, pYAP, and pAKT stains serve as assay controls. Bottom row: segmentation of the nuclei (empty circles) and dots (black spots). One out of two representative experiment is shown. Scale bar: 50 µm. (H) Quantification of the interactions between YAP/AKT and pYAP/pAKT in the nucleus (experiment shown in G). One symbol represents one image (YAP = 17, AKT = 16, YAP-AKT = 22, pYAP = 24, pAKT = 25, pYAP-pAKT = 23). (I) PLA showing YAP-AKT interaction units per nuclei of untreated (ctrl, n = 12) and APAP-treated cells (10 mM, for 2 hr, n = 20). APAP treatment induces nuclear interaction between YAP and AKT. Statistical test: unpaired two-tailed parametric t-test (p value = 0.0155). For A and C-F, actin served as loading control.
The online version of this article includes the following source data for figure 4: Source data 1. Source data for Western Blots shown in Figure 4A,C-F.

Figure 4 continued
were subjected to expression profiling and results were investigated regarding the presence of gene signatures specific for ROS, AKT, and Hippo/YAP activity (De Marco et al., 2017;Han et al., 2008;Wang et al., 2018). The experimental results showed that gene expression of ROS, AKT, and Hippo/ YAP target genes in livers was altered and prominently activated between 6 hr and 2 days after APAP treatment ( Figure 5A-C). To compare the temporal dynamic of the gene signatures, z-scores of genes in the signature per time point were summarized and normalized to the number of genes in the signature. As indicated by the results from the cell culture experiments, the expression index illustrated a specific order of signature activation starting with ROS (6-12 hr), followed by AKT (12 hr to 1 day) and Hippo/YAP signature genes (1-2 days) ( Figure 5D).
To confirm the findings from the gene expression analysis, we performed immunohistochemical staining of liver tissues isolated at five time points after APAP treatment. Staining for total YAP illustrated not only a general increase of YAP positivity in the cytoplasm of surviving hepatocytes, but also its prominent nuclear accumulation after 1-2 days ( Figure 5E, arrows). Moreover, nuclear AKT enrichment was already detectable 6 hr after APAP treatment. Interestingly, nuclear AKT accumulation in hepatocytes remained high for the rest of the experiment compared to untreated mice ( Figure 5E, arrows).
Because of a high mouse-to-mouse variability, we decided to objectify these results. For this, the immunohistochemical staining of YAP and AKT was quantitatively analyzed using machine learning (random forest) and image analysis methods ( Figure 5F). For this, 733 (for YAP) and 516 (for AKT) images were unbiasedly selected from all investigated animals and analyzed with an algorithm, which was trained to detect positively stained nuclei outside of the necrotic areas and tissue borders (Figure 5-figure supplement 1; for detailed information see Materials and methods). The quantification results confirmed the observed nuclear YAP positivity 1-2 days after APAP administration, while nuclear AKT was already detectable at 6 hr and maintained during the experiment ( Figure 5F). The shift from exclusive nuclear AKT to nuclear YAP was illustrated by direct comparison of the time points 6 hr and 2 days after APAP injection ( Figure 5G). The mechanistic connection between AKT activity and YAP induction in murine hepatocytes was confirmed in independent experiments. Here, hydrodynamic gene delivery of myristoylated AKT led to nuclear enrichment of YAP expression in hepatocytes and subsequent induction of typical YAP target genes (data not shown).
To sum up, we confirm the sequential order of APAP-induced events that cause YAP activation in mouse liver tissue.

Discussion
In this study, we aimed to investigate the different properties of the Hippo pathway effectors YAP and TAZ under in vitro and in vivo conditions that resembled the situation in DILI. For this, we developed a cell-based model system that allowed the systematic and comparative tracing of spatial YAP/TAZ dynamics upon stimulation with APAP, which is one of the most widely used analgesic with high potential of liver toxicity. By integrating time-resolved experimental data, computational modeling, image analysis tools as well as confirmatory in vivo results, we could draw several important conclusions on the role of YAP and TAZ under physiological conditions and in liver cells upon DILI.
Several mathematical modeling approaches have been previously applied to explain the dynamic behavior of the Hippo pathway under physiological or pathological conditions in different organisms. For example, computational modeling was used to explain how Hippo signaling cross-talks with other signaling pathways via protein interactions (Labibi et al., 2020;Romano et al., 2014). Other studies exclusively focused on biophysical/biochemical regulation of the Hippo pathway to investigate how YAP/TAZ activity relates to mechanical input information (Ege et al., 2018;Eroumé et al., 2021;Gou et al., 2018;Scott et al., 2021;Sun et al., 2016).
The goal of our approach was to determine whether nuclear phosphorylation of YAP plays a role in Hippo signaling, specifically in governing YAP localization. A model not including nuclear phosphorylation of YAP could not sufficiently explain the experimentally observed localization pattern and was excluded from further analyses. However, the alternative model proposed in this work, can explain the experimental data and therefore establishes a possible mechanism for YAP localization. However, more complex mechanisms regulating subcellular localization and dynamics of YAP and TAZ cannot be excluded. Therefore, our model partly explains how the Hippo pathway is organized, however, it does not aim to comprehensively describe it. Here, we decided to use PDE modeling, which allowed us to discover spatial aspects of signaling dynamics and to use live cell image data for parametrizing the Hippo signaling steady-state models. By doing so, we compared PDE model simulations of different model topologies with the experimentally acquired data for YAP as well as TAZ and selected a model topology, which best corresponded to data generated in vitro. The selected alternative Hippo pathway model is based on experimental observations and parametrized within the considered boundaries of experimental evidence (see Appendix 1). The obtained model parameters are results of the best model fit obtained, and allow reproducing the described model behavior. However, since the model fitting is based on steady-state data only and concentration observations are in arbitrary units, the parameter values are subject to structural non-identifiability. Nevertheless, this does not affect the qualitative conclusions drawn from the models for the investigated biological process.
The results of our PDE modeling approach strongly suggested that the nuclear phosphorylation reaction of YAP/TAZ was necessary and sufficient for the reproduction of their subcellular distribution patterns in vitro (i.e., the halo-like nuclear distribution of YAP and TAZ). This observation is of importance since the conventional view on the Hippo pathway topology represents sequential phosphorylation steps in the cytoplasm. Our PDE modeling-based finding would therefore extend the current understanding of Hippo pathway regulation since it adds the cell nucleus as a pivotal spatial compartment in the modulation and adjustment of Hippo/YAP/TAZ signaling. This conclusion is supported by our experimental findings illustrating that LATS1/2 are expressed in cell nuclei and that the interaction between phosphorylated LATS1/2 and YAP as well as between AKT and YAP is also detectable in this compartment. Importantly, previously published experimental and computational approaches support our findings. For example, nuclear LATS1/2 controls YAP phosphorylation in the context of NF2/Merlin deficiency (Li et al., 2014). In addition, the combination of quantitative photobleaching experiments and modeling illustrated that nuclear export is a crucial reaction for the subcellular YAP distribution (Ege et al., 2018). Combining these insights on nuclear export mechanisms with our findings on nuclear YAP phosphorylation strongly suggests that both processes closely cooperate in the efficient regulation of the Hippo pathway effectors. Lastly, our model does not explicitly exclude phosphorylation reaction in the cytoplasm; however, it suggests the necessity to reevaluate the canonical Hippo pathway scheme. The appeal for the reevaluation and extension of the canonical Hippo pathway, for example, by combining models on nuclear phosphorylation and nuclear transport, would further broaden our understanding on the dynamic spatial behavior of complex signaling pathways (Ege et al., 2018;Shreberk-Shaked and Oren, 2019).
The extension of the canonical Hippo pathway model with processes that control active nuclear YAP phosphorylation illustrates that the regulatory YAP/TAZ phosphorylation must be considered as continuum process (Shreberk-Shaked and Oren, 2019). For instance, the existence of nuclear phosphorylation possibilities by LATS-dependent and -independent (e.g., by AKT) mechanisms increases the complexity but also flexibility and redundancy of the signaling pathway (Ege et al., 2018;Gao et al., 2017;Low et al., 2014). These cellular aspects are part of a fine-tuned regulatory network that allows a rapid and adjustable proliferative response of YAP and TAZ under diverse physiological and pathological cellular conditions (e.g., upon induction of tissue regeneration).
To our knowledge, a direct mathematical comparison of YAP and TAZ regarding their spatial distribution has not been performed, yet. Using the PDE model, we comparatively investigated YAP and TAZ distribution and identified kinetic parameters, which discriminate between both species. Specifically, the observed localization differences (nuclear exclusion vs. nuclear localization) between YAP and TAZ might be caused by distinct phosphorylation rates of YAP and TAZ. Nuclear/cytoplasmic over time illustrate the timely order of ROS (max values between 6 and 12 hr), AKT (max values between 12 hr and 1 day) and YAP (max values between 1 and 2 days) target gene signatures after APAP treatment. (E) Mouse liver tissue sections stained for YAP (top tow) and total AKT (bottom row) under control condition and at 6 hr as well as 1, 2, 6, and 16 days after APAP treatment (300 mg/kg). Arrows indicate high YAP and AKT positivity in nuclei. Scale bar: 50 µm. (F) Automatic quantification of nuclei with YAP and AKT positivity in mouse liver tissues. Each dot represents the count of stained nuclei in one image (a tile, 1 mm 2 ) for YAP (n = 733) and AKT (n = 520). (G) Quantification of YAP and AKT nuclear positivity 6 hr and 2 days after APAP treatment. One dot represents the average count of positive nuclei in one animal (n = 4).
The online version of this article includes the following figure supplement(s) for figure 5:  shuttling dynamics for YAP and/or TAZ have been investigated previously (Ege et al., 2018;Plouffe et al., 2018). However, our model predicts that the phosphorylation efficiency for TAZ must be higher than for YAP to achieve the observed differential localization of YAP and TAZ.
By utilizing our image analysis approach for the detection of minor changes in the subcellular distribution of YAP and TAZ, we showed a dynamic shuttling response of YAP (and to a much lesser extent of TAZ) upon APAP treatment. Although, APAP-dependent activation of YAP was described in the literature (Poudel et al., 2021), the underlying molecular mechanisms were poorly understood. Our results demonstrated that APAP-induced YAP shuttling is not the result of an unspecific cellular response but is due to a sequential order of molecular events that include ROS induction (Ghallab et al., 2016;Shuhendler et al., 2014), followed by AKT (in)activation (Koundouros NPoulogiannis, 2018) and nuclear YAP (de)phosphorylation. Time-dependent APAP-induced ROS activation was recently supported by a study by Ghallab et al., which showed transient induction of oxidative stress in mice 8 hr after APAP overdose with a peak at 2 hr (Ghallab et al., 2022). Interestingly, we observed a biphasic YAP regulation: YAP phosphorylation/inactivation at early time points (up to 3 hr) and its dephosphorylated/activation at later time points (>24 hr). Especially, the functionally important activation of YAP, which is associated with the induction of cell proliferation, is well reflected by the temporal activation of YAP in murine tissues upon APAP administration after 24-48 hr.
Importantly, several cellular mechanisms might simultaneously contribute to the APAP-dependent activation of YAP. For example, the direct impact of APAP on MST1/2 kinases, which are central Hippo pathway constituents, represents one additional process that may contribute to the observed effects on YAP (Fan et al., 2016). Moreover, the serine-/threonine-kinase AKT could control the Hippo/YAP axis at different levels. For example, a physical interaction of AKT and YAP in the nucleus (shown here) or the cytoplasm might directly control YAP phosphorylation (Basu et al., 2003). Alternatively, AKT can phosphorylate MST2 and therefore indirectly contribute to YAP activity (Romano et al., 2014). Thus, APAP probably controls YAP in a multimodal manner to achieve a cellular response and it is therefore likely that this multimodal process is part of cell-protective mechanism that counteracts the massive loss of hepatocytes upon APAP-induced DILI (Holland et al., 2022;Shuhendler et al., 2014).
In our manuscript, we show that APAP administration induces a chain of molecular events through the APAP/ROS/AKT axis, which is leading to YAP inactivation (early) followed by YAP activation (late). The late induction of YAP activity might be considered as cell-protective cellular response upon longterm tissue damage; however, YAP also acts as a potent oncogene in different cell types, including hepatocytes, and its overexpression is associated with tumor initiation (Dong et al., 2007;Weiler et al., 2017). Thus, it is tempting to speculate that APAP overdose may not only cause DILI but also increases the risk for liver tumor development. Although controversially discussed in the literature, a recent evaluation of 139 published epidemiologic studies strongly argues against a relationship between APAP exposure and cancer (Weinstein et al., 2021). Actually, APAP has been discussed as therapeutic option for patients with YAP-induced cancer (Poudel et al., 2021). Although, our and other studies clearly demonstrate a mechanistic link between APAP uptake and Hippo/YAP pathway activity, these partly controversial findings illustrate the necessity to decipher this connection under distinct disease conditions in the future.

Western immunoblotting
Western immunoblotting was performed as previously described (Weiler et al., 2017). In brief, for total protein fraction preparation, cultured cells were harvested using 1× Cell Lysis Buffer (Cell Signaling Technology, Danvers, USA) supplemented with 1× Protease Inhibitor Mix G (Serva, Heidelberg, Germany) and 1× PhosSTOP (Roche Diagnostics Deutschland GmbH, Mannheim, Germany). For protein fractionation, the NE-PER Nuclear and Cytoplasmic Extraction kit was used (Thermo Scientific, Rockford, USA). The following antibodies were used for western immunoblotting experiments: Western blot detection and quantification were performed using the Odyssey-CLx Infrared Imaging system with the ImageStudio software (LI-COR Biosciences, Bad Homburg, Germany). Phosphorylated protein band intensity was measured and normalized to total protein concentrations. Equal amount of protein was loaded for each line, as measured by the Bradford reagent (Millipore Sigma, Saint Louis, USA). Raw unedited image files and uncropped blots are available as source data files of this manuscript.
RNA isolation, reverse transcription, and quantitative real-time PCR (qPCR) RNA was isolated using Extractme kit (Blirt, Gdańsk, Poland) according to the manufacturer's protocol. Reverse transcription was performed using the RevertAid kit (Thermo Scientific). Semiquantitative real-time PCR reactions were set up using the primaQuant 2× qPCR-SYBR-Green-Mastermix (Steinbrenner Laborsysteme, Wiesenbach, Germany) and analyzed with the QuantStudio 3 real-time PCR system (Applied Biosystems, Thermo Fisher Scientific, Singapore). The following cycling conditions were applied: 95°C for 15 min, followed by 40 cycles of 95°C for 15 s, and 60°C for 60 s. Product specificity was confirmed by melting curve analysis (95°C for 15 s, 60°C for 30 s, 60-95°C with 0.5°C/s).
The mRNA levels were normalized to glyceraldehyde-3-phosphate dehydrogenase (GAPDH), 60S ribosomal protein L41 (RPL41), serine/arginine-rich splicing factor (SRSF4), and β2-microglobulin (B2M). The following primers for human cDNAs were used: YAP-for: For the analysis of gene expression in tissue samples, a panel of housekeeping genes was analyzed using the geNorm algorithm to find the most stable reference gene (Vandesompele et al., 2002).

In situ PLA
The DuoLink in situ PLA was performed according to the manufacturer's instructions (Sigma-Aldrich). Briefly, cells were seeded on glass coverslips and prior to APAP administration cells were grown under FCS-free conditions for 1 day, then incubated with APAP (10 mM) or PBS. After treatment, cells were washed three times with 2 mM MgCl 2 in PBS and fixed with 4% paraformaldehyde for 10 min at room temperature. Fixed cells were washed four times with PBS for 5 min, permeabilized with 0.2% Triton X-100 in PBS for 5 min at room temperature, and washed again twice with PBS for 5 min. Subsequently, cells were blocked with Blocking solution (Sigma-Aldrich) for 30 min at room temperature and incubated with the following primary antibodies diluted in Antibody Diluent (Sigma-Aldrich) overnight at 4°C: anti-YAP (1:25, Santa Cruz Biotechnology, sc-271134, RRID:AB_10612397), anti-AKT (1:200, Cell Signaling Technology, #9272, RRID:AB_329827), anti-pYAP (1:200, Cell Signaling Technology, #4911, RRID:AB_2218913), anti-pAKT (1:200, Cell Signaling Technology, #4051, AB_331158), and anti-pLATS1/2 (1:200, Cell Signaling Technology, #8654, RRID:AB_10971635). Subsequently, cells were washed twice with Wash Buffer A (Sigma-Aldrich) and incubated with prediluted rabbit PLUS and mouse MINUS probes (Sigma-Aldrich) in Antibody Diluent for 1 hr at 37°C. After incubation, cells were washed twice with Wash Buffer A and then incubated with ligation solution (Sigma-Aldrich) for 30 min at 37°C. After ligation, samples were again washed twice with Wash Buffer A for 2 min at room temperature and incubated with amplification solution, and detection reagent Orange (Sigma-Aldrich) for 100 min at 37 °C. Finally, cells were washed twice with Wash Buffer B (Sigma-Aldrich) for 10 min at room temperature, once with 0.01× Wash Buffer B for 1 min at room temperature and coverslips were mounted on the slide with DAPI Fluoromount-G mounting medium (SouthernBiotech, Birmingham, USA).
Fluorescence images were captured using an inverted Nikon Ti2 microscope with Nikon S Plan Fluor ELWD ×40 NA 0.60 objective in a widefield fluorescence mode using Lumencor Sola SE II lamp. Images were captured in DAPI and TRITC channels (460 and 580 nm) with Nikon DS-Qi2 monochrome camera (image size 2404 × 2404 px, pixel resolution of 0.18 µm/px).

ROS activity measurement
Measurements of the intracellular ROS levels were performed using the DCFDA/H2DCFDA cellular ROS assay kit (Abcam, Amsterdam, Netherlands) according to the manufacturer's protocol. In brief, cells were seeded on white clear-bottom 96-well plate (Corning, Corning, USA) and incubated overnight. Cells were treated with 10 mM APAP, 2 mM H 2 O 2 , or 300 µM TBHP (H 2 O 2 and TBHP served as positive controls for ROS induction) for 6 hr in FCS-free cell culture medium. Fluorescence was measured using a microplate reader (FluoStar Omega, BMG Labtech GmbH, Ortenberg, Germany). Buffer solution without cells served as a background control.
In vivo experiments and sample analyses Housing and treatment of mice and induction of acute liver injury by acetaminophen Male C57BL/6N mice (8-to 10-week-old) were bought from Janvier Labs (Janvier Labs, Le Genest-Saint-Isle, France). Animals were housed under 12 hr light/dark cycles at controlled ambient temperature of 25°C with free access to water and were fed ad libitum with a standard diet (Ssniff, Soest, Germany) before starting the experiments. Induction of acute liver injury with APAP was done as previously described (Schneider et al., 2021). Briefly, the mice were fasted overnight, then challenged with a dose of 300 mg/kg APAP intraperitoneally. APAP was dissolved in warm PBS with an application volume of 30 ml/kg. Control group was treated with PBS only. The mice were fed ad libitum after APAP administration. All animals were included for further analyses. All experiments were approved by the local animal welfare committee (LANUV, North Rhine-Westphalia, Germany, application number: 84-02.04.2016.A279). Hydrodynamic gene delivery experiments in mice were performed as recently described (Luiken et al., 2020).

Expression profiling and bioinformatics
RNA isolation and gene array analysis were performed as published before (Campos et al., 2020) and bioinformatic analysis was done as described (Holland et al., 2022). The gene expression data are available under ArrayExpress accession number GSE167032. The expression data were applied to three known signatures that are informative for ROS activity (Han et al., 2008), AKT activity (De Marco et al., 2017), and YAP/TAZ activity (Wang et al., 2018). Due to the high number of AKT signature genes, only genes whose response to APAP treatment was larger than fold change of 2 (compared to control animals) were considered. The expression data were z-score normalized and clustered with seborn clustermap python module (v0.11.0). The signature score was obtained by summarizing z-scored expression values at the given time point if the z-scored value was greater than 0.5. The summarized expression values were normalized to the number of genes in the signature.

Analysis of the live cell images
The confocal images of the living cells were analyzed in a high-throughput manner using ImageJ (v1.53f51) platform . Images were first manually selected with respect to the quality criteria: images with insufficient sharpness or with artifacts were discarded from the analysis. The image processing pipeline was based on Weka segmentation (v3.3.1) (Arganda-Carreras et al., 2017) of foreground and background areas, and subsequent thresholding, object detection, and counting. After object detection and counting, the respective masks were overlaid on the initial images to acquire YAP and TAZ intensity values for nuclei and cytoplasm. Mean pixel intensity from nuclear areas of the cells was divided by the mean pixel intensity of cytoplasmic regions, thus obtaining an NCR.

Analysis of PLA images
PLA slides were analyzed using ImageJ (v1.53f51). First, nuclei and dots were classified using the Weka segmentation algorithm (v3.3.1) (Arganda-Carreras et al., 2017), thresholded, and counted. The pseudo-cytoplasmic (ring-shaped) area was created using the ImageJ's binary mask option dilate for 30 iterations on nuclear masks to obtain nuclear and cytoplasmic area for a comparative analysis. The thresholded nuclei or cytoplasmic masks were overlaid with detected dots to obtain the information on the subcellular localization of the protein interaction.
First, digital images of the tissue sections were divided in tiles (1 mm 2 , 400 × 400 pixels), which were binned (to reduce processing time and storage load) using ASAP and python modules. Some images were excluded due to staining and/or scanning artifacts. Tiles, which displayed tissue for at least 50% of their area, were kept for further processing. Machine-learning model, which was based on a random forest algorithm, was trained on a selected set of training tiles for YAP or AKT using Ilastik software. The algorithm was trained to detect positively stained nuclei, excluding necrotic areas and staining artifacts. Ilastik software was further used to export probability maps, which were possessed with ImageJ. In ImageJ, probability maps were thresholded and positive nuclei were counted and normalized to the area of the tile, which was occupied by the tissue.

PDE modeling
Spatial modeling aimed at describing the variations in space of fluorescently labeled YAP and TAZ distribution within living Hep3B cells. The mathematical model is based on a system of PDEs. PDE modeling and parameter estimation were performed with the Spatial Model Editor (SME) software (v1.2.1) and sme-contrib (v0.0.14) python module (https://spatial-model-editor.github.io/). SME is graphical user interface-based model editing and simulation software compatible with systems biology markup language (SBML) standards. Models were simulated using simple Forward Time Centered Space (FTCS) solver.
In the demonstrated models, all reactions were defined by first-order kinetics (for PDEs see Appendix 1). For the parameter estimation, the model behavior was evaluated at steady state, that is the model was simulated until the concentration distribution did not change over time. Parameter estimation was performed using the particle swarm algorithm (20 particles, 200 iterations) to minimize a cost function consisting of the weighted sum of two terms: the squared per pixel differences between model and the target image, and the sum of squares of species concentration rates of change. In total 200 fitted parametrizations for each tested model were generated. The parameter space for the optimization algorithm was defined based on published data (Appendix 1). The parameter ranges of our canonical and alternative models were selected in a way to describe biologically meaningful value ranges and to avoid numeric instability during PDE simulations of the canonical model. The mathematical models were uploaded to the BioModels repository under the model identifier number MODEL2202080001. PDE model equations and parameters can be found in Supplementary Information (Appendix 1-Tables 1-6; Ege et al., 2018;Jack et al., 1990). Tables showing all fitted parameters for YAP and TAZ are provided (Supplementary file 1, Supplementary file 2).

Data availability
The mathematical models generated and analyzed during the current study are available in the BioModels repository (https://www.ebi.ac.uk/biomodels/) under accession number MODEL2202080001. The PDE simulation software is available over Github (https://spatial-model-editor.github.io/). The gene expression data analyzed during the current study are available in the Gene Expression Omnibus (GEO) repository under accession number GSE167032.
The following dataset was generated: